Charge Order in a Two-Dimensional Kondo Lattice Model 



m 
o 



Takahiro Misawa, Junki Yoshitakc, and Yukitoshi Motomc 
Department of Applied Physics, University of Tokyo, 
7-3-1 Hongo, Bunkyo-ku, Tokyo, 113-8656, Japan 
(Dated: March 14, 2013) 

The possibility of charge order is theoretically examined for the Kondo lattice model in two 
dimensions, which does not include bare repulsive interactions. Using two complementary numerical 
methods, we find that charge order appears at quarter filling in an intermediate Kondo coupling 
region. The charge ordered ground state is an insulator exhibiting an antiferromagnetic order at 
charge-poor sites, while the paramagnetic charge-ordered state at finite temperatures is metallic 
with pseudogap behavior. We confirm that the stability of charge order is closely related with the 
local Kondo-singlet formation at charge- rich sites. Our results settle the controversy on charge order 
in the Kondo lattice model in realistic spatial dimensions. 

PACS numbers: 71.10.Fd, 71.27.+a, 75.25.Dk, 71.30.+h 
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Charge order (CO), a periodic spatial modulation of 
electron density which breaks the lattice translational 
symmetry, is a fundamental electronic state of matter. 
Since it was argued for the Verwey transition in mag- 
netite [H, CO has been ubiquitously found in a broad 
range of correlated electron systems, such as organic 
conductors [H, transition metal compounds 0], and f- 
electron compounds 0, [Bj] . One of the obvious origins of 
CO is the intersite Coulomb interaction, which prevents 
electrons from coming close to each other. For instance, 
a phase transition from paramagnetic (PM) metal to CO 
insulator takes place in an extended Hubbard model at 
quarter filling while increasing nearest-neighbor repulsive 
interaction [ff. 

An alternative, interesting possibility of CO was pro- 
posed for a Kondo lattice model which does not in- 
clude any bare repulsive interaction between electrons 
Eq. Q] 0. The CO instability was deduced from a 
perturbation expansion in the limit of large onsitc Kondo 
coupling J. The second-order perturbation in terms of 
the electron hopping t gives rise to an effective repulsion 
between neighboring electrons oc t 2 /J, whose competi- 
tion with the kinetic energy oc t was expected to cause 
CO at quarter filling in an intermediate t/J regime. After 
the interesting proposal, the ground state of the Kondo 
lattice model in one dimension was studied numerically, 
but such CO instability was not found Recently, 
however, the problem was reexamined in infinite dimen- 
sions by employing the dynamical mean-field theory, and 
the instability toward CO was found near quarter fill- 
ing [ll|, Thus, the conclusions arc contradicting 
between one and infinite dimensions, and it is unclear 
whether CO is realized in more realistic two and three 
dimensions. This is a fundamental issue related with the 
unique mechanism of CO in the Kondo lattice model, not 
seen in Hubbard-type models, and potentially relevant to 
understand CO phenomena in /-electron compounds. 

In this Letter, we explore CO in the quarter-filled 
Kondo lattice model on a square lattice by using two com- 



plementary methods: variational Monte Carlo (VMC) 
method [1 31 ] for the ground state and cellular dynamical 
mean-field theory (CDMFT) [3l for finite temperatures 
(T). In the VMC simulation, the optimization of large 
number of variational parameters and careful system-size 
analysis allow us to obtain the ground state by seriously 
taking account of both spatial and quantum correlations. 
On the other hand, in CDMFT, we employ a continuous- 
time Monte Carlo method [IH for the impurity solver, 
which enables us to obtain finite- T properties with fully 
including dynamical correlations as well as spatial corre- 
lations within the cluster. The results by the two meth- 
ods provide convincing evidence for the appearance of 
CO. Furthermore, we show that the CO ground state 
is insulating and accompanied by an antiferromagnetic 
(AF) order. At finite T, the charge gap closes as AF 
order disappears, while a pseudogap remains in the PM 
CO phase. We confirm that CO is stabilized by the lo- 
cal Kondo-singlet formation by studying the effect of an 
intersite AF exchange interaction. Our results suggest 
that the Kondo lattice model accommodates CO in more 
than one dimension. 

We consider the Kondo lattice model on a square lat- 
tice, whose Hamiltonian is given by 

H = "* E (4,<*r + h-C) + S > ■ S i ' W 

{*,]),<* i 

where c\ a (c^) is a creation (annihilation) operator of 
a conduction electron with spin a at ith site. The 
spin operator of conduction (localized) electron at ith 
site is defined by S- = \ Aa ™' ™ 1 ( S i = 

sE^'/l^/w'i fL and fi<r are creation and an- 
nihilation operators for the localized electrons, respec- 
tively), where <j aa i is a vector of the Pauli matrices. The 
first term describes the hopping of conduction electrons 
between the nearest-neighbor sites on the square 

lattice, and the second term represents the onsite AF 
Kondo coupling between conduction and localized spins 
(J > 0). In the following, we set the lattice constant 
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FIG. 1: (color online), (a) Finite-T phase diagram for the 
model (O at quarter filling obtained by CDMFT. Tco and 
Taf denote the critical temperatures for charge ordering and 
antiferromagnetic ordering, respectively, (b) T dependencies 
of the local charge density in the A and B sublattice, ua 
and ria, at J/t — 3. (c) T dependencies of the magnitude of 
local magnetic moment of conduction electrons in the A and 
B sublattice, tua and Mb, at J/t = 3. The lines connecting 
the data points are the guides for eyes. (d)(e) Sublattice for 
the two- and four-site cluster calculations. See the text for 
details, (f) Schematic picture of the CO phase with AF order 
at the lowest T. The size of circles represent the magnitude 
of local density and the arrows show the magnetic moments. 



a = 1 and focus on the quarter filling case, namely, at 
n c = ^£ 1<T ( C L C ^) = 1/2 (N s = L x L is the system 
size). 

We study the possibility of CO in the model ([T]) by 
two complementary methods, VMC and CDMFT. In the 
VMC calculations, we study the ground state proper- 
ties by employing a generalized BCS wave function with 
the quantum number projection and the Gutzwiller and 
Jastrow factors; = Vg^jC 3 ]^^). Here, C s is 
the spin projection operator to the total spin S sub- 
space; Vg and Vi are the Gutzwiller and Jastrow fac- 
tors, respectively [13j. The spin projection is per- 
formed onto the S = singlet subspace, except for 
the case explicitly mentioned. The Gutzwiller factor 
penalizes the double occupation of electrons by Vg = 
CX P(- Ei.A^n^n^) where n x ia = \\ a X i<T (A = c,/); wc 
take g{ — oo for localized / electrons. The Jastrow fac- 
tor is introduced only for conduction electrons as Vj = 

exp(-|Ei,j u ij<«i)) where n i = T,a n i<x- The one- 
body part Impair) is the generalized pairing wave function 

defined as |0 pair ) = (Ex^ cJ E£=i /^Aj t ^)^/ 2 |0), 
where N e is the number of electrons (including / elec- 
trons). In this study, we restrict the variational pa- 
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structure. All the variational parameters are simultane- 
ously optimized by using the stochastic reconfiguration 
method OUB]. Our variational wave function \tp) can 
flexibly describe CO, AF, and PM state on equal footing. 
The calculations are done up to 12 x 12 sites; to reduce the 
finite-size effects, wc choose appropriate boundary condi- 
tions for each system size so as to satisfy the closed-shell 
condition in the noninteracting case J = 0. In addi- 
tion to the ground state by VMC, we calculated finite-T 
properties by CDMFT for two-site (1x2) and four-site 
(2 x 2) clusters [see Figs. [TJd) andQJc)]. For the two- 
site case, we effectively consider a 2 x 2 sublattice order 
by configuring the Green's function with flipped spins 
so as to accomodate the ordering pattern found in the 
four-site cluster [Fig. Hfd)] . To solve the effective impu- 
rity problems for the clusters, we employ the continuous- 
time Monte Carlo technique [ll|. Typically, we repeat 
the CDMFT loops for 50 times for convergence, and per- 
formed 10 6 -10 7 Monte Carlo measurements per loop. 

Let us first discuss finite-T properties. Figure Ufa) 
shows the finite-T phase diagram obtained by the 
CDMFT calculations. The result shows that the model 
in Eq. ((TJ) exhibits CO at quarter filling in the interme- 
diate coupling region for 1.5 < J/t < 3.5. CO is a two- 
sublattice checkerboard type, as shown in Fig. [ljf). T 
dependence of the charge disproportionation at J/t = 3 
is shown in Fig.JIJb). The continuous development of the 
order parameter tia — «b as lowering T suggests that this 
CO transition is of second order. As shown in Fig.QJa), 
the critical temperature Tco grows rapidly as increasing 
J/t, whereas it suddenly disappears for J/t > 3.5 in the 
T range calculated. The behavior of Tco is very similar 
to that obtained in infinite dimensions (corresponding to 
a single-site cluster calculation) [ll|. Interestingly, Tco 
does not largely depend on the cluster sizes in the present 
calculations. 

In the CO phase, an AF order appears in the low T 
region. In Fig. [He), wc show T dependences of the local 
magnetic moments of conduction electrons at the charge- 
rich site, mA, and charge-poor site, roe- At the charge- 
poor site, a spin polarization develops continuously as 
lowering T, while no magnetic moment is seen at the 
charge-rich site. The AF order is a four-sublattice one, as 
shown in Fig.JIJf); collinear Neel type magnetic moments 
appear at the charge-poor sites, while the charge-rich 
sites remain nonmagnetic. Strictly speaking, the crit- 
ical temperature Taf should be zero, as the model has 
SU(2) symmetry in two dimensions [17] : the finite Taf is 
an artifact of CDMFT. The result, however, reflects the 
development of AF correlations toward T = 0. Hence, 
the CDMFT results suggest that the ground state is the 
CO state with AF long-range order for 1.5 < J/t < 3.5. 

Now, we investigate the ground state by the VMC cal- 
culations. Figures HJa) and [U[b) show J/t dependences 
of the peaks of the charge structure factor and spin struc- 
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FIG. 3: (color online). Chemical potential (i as a function 
of the electron density n c at J/t — 3. Gray lines in (a) are 
the guides for eyes. The jump of \x at quarter filling n c = 0.5 
corresponds to the charge gap A c in the CO state. CDMFT 
results are for the two-site cluster. The inset of (a) shows 
the momentum distribution function n(fc) at J/t = 3 and 
n c = 0.5 obtained by VMC. 



FIG. 2: (color online). J/t dependences of the peak values of 
(a) charge structure factor and (b) spin structure factor for 
system sizes L = 6, 8, 10, and 12 calculated by VMC. System 
size extrapolations of the data are shown in (c) and (d). The 
extrapolated values in the bulk limit L —¥ 00 are plotted in 
(a) and (b). 

ture factor, which arc defined as 

s i,j 

s ^=w ^ {{s * + s * } ■ + s ')> eM * <_r,) > ( 3 ) 

respectively. Here, r, is the positional vector for the ith 
site. We find that N(q) has a single peak at q pcak = 
(ir, 7r) for all the system sizes in the entire region of J/t. 
On the other hand, S(q) shows two equivalent peaks at 
9pcak = C 71 ") 0) an d (0, 7r) simultaneously, whereas it shows 
a single peak at (tt, 0) or (0, ir) for small system sizes 
and for J/t < 2 [shown by crosses in Fig. [2fb)]; we plot 
£(<? P eak) = {5(tt, 0) + 5(0, tt)}/2. The data well scale 
with 1/L as shown in Figs. [U^c) and[3{d); the extrapo- 
lation to L — >• cxd indicates that CO accompanied by AF 
order appears in the bulk limit for J/t > 2. The size ex- 
trapolation of S(q) is not shown for J/t < 2 because the 
magnetic ordering patterns depend on system sizes. To 
clarify the magnetic order in this region, it is necessary to 
extend the sublattice for variational parameters; this is 
left for a future study. The VMC data indicate that the 
two-sublattice CO continuously appears for J/t > 1.5, 
accompanied by the four-sublattice AF order, as illus- 
trated in Fig. Qif). 

In the large J/t region where Tco disappears in 
CDMFT, we find that the VMC solution shows an in- 
stability toward a ferromagnetic state. Using the spin 
projection, we calculate the ground state energy as a 



function of the total spin 5, E(S). We found that E(S) 
has local minima at 5 = and 5 = iV s /4 for large J/t, 
which correspond to the singlet and ferromagnetic states, 
respectively. In the ferromagnetic state, conduction elec- 
trons form local singlets with localized spins at each site, 
and the rest unpaired localized moments are fully polar- 
ized. The result shows that E(S = 0) < E(S = N s /A) for 
J/t < 4, but the energy difference becomes smaller as in- 
creasing J/t- E(S = 0) ~ E(S = N s /4) for J/t ~ 5 [ij]. 
This implies that the system undergoes a transition from 
the CO+AF state to the ferromagnetic state at J/t ~ 5. 

Therefore, our complementary study by CDMFT and 
VMC clearly indicates that the CO+AF ground state 
appears in the intermediate J/t region for 1.5 < J/t < 
3.5 in two dimensions. The result is in sharp contrast to 
the recent one obtained in infinite dimensions, in which 
the CO state appears only at finite T and is replaced by a 
charge-uniform ferromagnetic state at low T in the entire 
range of J / 1 [l2| ■ 

Our results suggest that CO is also stabilized in three 
dimensions because spatial fluctuations are suppressed 
as the spatial dimension becomes larger. It is, however, 
nontrivial what type of magnetic order appears in the 
(7r,7r,7r)-CO state, as the magnetic moments at charge- 
poor sites suffer from geometrical frustration of the face- 
centered-cubic type. The interesting problem will be dis- 
cussed elsewhere (l9j . 

Next, to clarify the electronic structure of the obtained 
CO state, we calculate the chemical potential fj, as a 
function of the electron density n c . In the VMC cal- 
culations, we calculate /i by using the relation n(n c ) = 
{E{n\) - E{n c 2 )}/{n<i - n|), where E{n c ) is the total 
energy at n c , and h c = (n\ + n|)/2. To reduce the 
finite-size effects, we choose n\ and n% that satisfy the 
closed-shell conditions in the noninteracting case. The 
results at J/t = 3 are shown in Fig. EJa). We find that 
fi shows a clear jump at quarter filling n c = 0.5. This 
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indicates that the CO+AF ground state is an insulator. 
The charge gap is estimated by the jump, A c ~ OAt at 
J/t = 3. The insulating nature is also observed in the 
momentum distribution function, calculated by 



s i,j,<r 



(4) 



In the metallic state, n(k) shows a discontinuity at the 
Fermi surface, while it does not show such singularity in 
the insulating state. The result for J/t = 3 at quarter 
filling is shown in the inset of Fig. [3fa) . n(k) shows a 
smooth k dependence, which supports that the system is 
insulating in the CO+AF ground state. 

Correspondingly, in the CDMFT results, /j, shows 
plateau-like behavior in the lowest-T CO+AF phase, as 
shown in Fig. [3]Ja). Here, the CDMFT calculations are 
performed in the grand canonical ensemble, and hence, 
n c is calculated for a given fi. As shown in Fig.^b), the 
plateau-like behavior is smeared out as raising T in the 
intermediate-T PM CO phase. The result suggests that 
the charge gap closes according to the disappearance of 
AF order. The PM CO phase at finite T is metallic, while 
pseudogap behavior remains as kink-like behavior with 
an inflection point in the fi-n data shown in Fig. [3{b). 
The conclusion is apparently in contradiction with the 
recent result obtained in infinite dimensions, where an 
insulating PM CO state is predicted [12 1. 

Finally, let us examine the origin of CO. Originally, 
the possibility of CO in the Kondo lattice model was sug- 
gested by the second-order perturbation theory in terms 
of t/J [7|. In this argument, the local Kondo spin-singlet 
formed in the limit of J/t — > oo plays a key role in the 
stability of CO. This picture was also examined by an 
energetic argument considering the Kondo effect [13] ■ To 
confirm the importance of local spin-singlet formation, 
we numerically study the effect of a perturbation which 
destroys the local singlets. Here, we consider a perturba- 
tion to the model ([I]) by adding the AF exchange inter- 
action between the localized spins, which is given by 



n AF = Jaf s ( 

<i,i) 



4 



(5) 



This term is expected to induce a Neel type AF or- 
der, which competes with the local spin-singlet forma- 
tion. The VMC results while changing Jaf A are shown 
in Fig. [4j As expected, the system exhibits a phase tran- 
sition from the CO+AF state to Neel AF ordered state 
without CO at Jaf A — 0.42. The transition point is es- 
timated by the comparison of the ground state energy for 
two states. The transition is of first order; the peaks of 
S(q) suddenly change from (tt, 0) and (0, 7r) to (tt, tt), and 
simultaneously, the peak of N(q) at (tt, 7t) is suppressed, 
as shown in Fig. [Ha). At the same time, as shown in 
Fig. 01b), the local spin correlation (Sf ■ S{), which is a 
measure of the local singlet formation, is also suddenly 
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FIG. 4: (color online), (a) Jaf dependences of the peak values 
of spin and charge structure factors at J/t = 3. (b) Local spin 
correlations at the charge-rich site (A) and charge-poor site 
(B). The data are obtained for iV s = 12 x 12 by VMC. 



suppressed. In particular, the local singlet at the charge- 
rich site is strongly suppressed at the transition. The re- 
sults indicate that the local Kondo-singlet formation at 
the charge-rich sites plays an important role to stabilize 
CO, in consistent with the perturbation argument. 

To conclude, using two complementary numerical 
methods, VMC at T = and CDMFT at finite T, we 
have provided convincing evidences for CO in the Kondo 
lattice model in two dimensions. CO appears in the in- 
termediate J/t region, accompanied by a collinear AF 
order at the charge-poor sites (the charge-rich sites are 
nonmagnetic). We also found that the CO+AF state is 
insulating, whereas the PM CO state at finite T is metal- 
lic while showing pseudogap behavior. We also confirmed 
that the local Kondo-singlet formation at the nonmag- 
netic sites plays an important role in stabilizing CO. Our 
results settle the controversy on the existence of CO in 
the Kondo lattice model in realistic dimensions. 

The present CO is of kinetic origin, and qualitatively 
different from the usual ones driven by bare intersite 
Coulomb repulsions. In the realistic compounds, how- 
ever, the intersite Coulomb repulsions are also present, 
which may cooperatively work with the Kondo coupling 
on the tendency toward CO. Nevertheless, it will be pos- 
sible to examine which mechanism is dominant by, e.g., 
applying magnetic field and pressure that affect the band- 
width, local Kondo coupling, and intersite interactions in 
a different way. 
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